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ABSTRACT 


The  bivariate  normal  distribution  function  may  be  expressed  as  the 


product  of  a  marginal  normal  distribution  times  a  conditional  distribution. 
By  approximating  this  conditional  distribution,  we  obtain  a  simple  method 
for  approximating  bivariate  normal  probabilities.  When  the  correlation 
falls  in  the  interval  [-.5,  .5],  the  maximum  absolute  error  in  our  ap¬ 
proximation  is  always  less  than  .0008.  The  conditional  distribution  that 
we  approximate  is  referred  to  as  a  ’normal  conditioned  on  a  truncated 
normal’  distribution  and  is  related  to  screening  and  selection  problems. 


1  ■  y  '  ^  .  >A..«  _ 


S'?' -  - £ 


i  ^v„w  r  u^'-iVAv  'n- 


INTRODUCTION 

The  bivariate  normal  distribution  is  frequently  employed  as  a  model 
for  screening  procedures,  where  a  product  is  either  accepted  or  rejected 
based  on  some  secondary  measurement  that  is  correlated  with  the  performance 
characteristic  (see,  e.g.,  Owen,  Mclntire  and  Seymour,  1975).  When  screen¬ 
ing  is  utilized,  the  proportion  of  all  items  which  are  rejected  even  though 
they  are  good  and  the  proportion  of  accepted  parts  which  are  in  fact  de¬ 
fective  are  of  interest.  The  purpose  of  this  paper  is  to  provide  simple 
formulae  for  approximating  these  probabilities.  The  results  that  follow 
should  also  be  of  interest  in  other  settings  where  bivariate  normal  prob- 
abilities-are  required  but  the  potential  user  does  not  have  ready  access 
to  a  computer. 

Tables  of  the  bivariate  normal  distribution  are  like  those  of  the  uni¬ 
variate  normal  in  that  the  variates  are  standardized,  i.e.,  each  variate 
has  its  mean  subtracted  and  the  difference  is  divided  by  the  standard  devia¬ 
tion.  Hence,  we  start  with  X  and  Y  having  a  joint  (standardized)  bivariate 
normal  distribution  with  correlation  p.  We  propose  a  simple  approximation 
for  the  cumulative  distribution  function  (cdf)  F(*,  ♦;  o)  of  (X,  Y).  For 
any  constants  h  and  k,  the  probability  F(h,  k;  p)  may  be  factored  into  a 
marginal  probability  times  a  conditional  probability,  i.e., 

F(h,  k;  p)  =  Pr[X  <  h]  •  Pr[Y  <  k  J  X  <  h].  (1) 

We  propose  approximating  the  conditional  distribution  of  Y,  given  X  <  h,  by 
a  normal  distribution.  In  this  way,  F(h,  k;  p)  can  be  approximated  as 
the  product  of  two  univariate  normal  probabilities. 

Mallows  (1959)  obtained  an  approximation  for  F(h,  k;  p)  which  requires 
simple  computations  and  evaluation  of  univariate  normal  probabilities  and 


’fSRKK'. 


I 


r-g^gggsgT 


percentiles.  For  the  values  of  h,  k  and  p  that  Mallows  examined,  the  maxi¬ 
mum  error  was  .0076.  The  approximation  we  present  is  simpler  than  Mallow's 
formulae  and  is  more  accurate  unless  jp|  exceeds  .8.  El  Lozy  (1982)  provides 
a  simple  way  to  compute  probabilities  of  the  special  form  F(h,  ph;  p).  El 
Lozy's  method  is  most  accurate  for  p  close  to  1 .  In  contrast,  our  approxima¬ 
tion  has  its  greatest  accuracy  when  jpj  is  less  than  .5. 

Before  we  present  the  approximation,  we  investigate  the  conditional  dis¬ 
tribution  of  Y,  given  X  <  h.  This  distribution,  which  we  label  the  normal 
conditioned  on  a  truncated  normal,  is  of  interest  in  itself  since  it  can 
represent  the  distribution  of  the  performance  variable  in  the  accepted  popula¬ 
tion  after  screening  with  respect  to  an  upper  specification  limit. 

NORMAL  CONDITIONED  ON  A  TRUNCATED  NORMAL 

It  is  well-known  that  the  distribution  of  Y  conditioned  on  X  =  h  is 

2 

normal  with  mean  ph  and  variance  1  -  p  .  That  is,  Y  is  said  to  have  a 
conditional  normal  distribution.  However,  the  distribution  of  Y  condi¬ 
tional  on  X  <  h  is  not  normal  (unless  p  =  0),  but  has  a  probability  density 
function  (pdf)  given  by 

G' (y)G[(h-py)/(l-p^)'i]/G(h) ,  ^2) 

where  G ‘ ( - )  and  G(*)  denote  the  standard  normal  pdf  and  cdf,  respectively. 

;s 
(3) 


Let  Y^  denote  the  random  variable  with  pdf  (2).  The  mean  of  Y^  is 


u  =  -pG*(h)/G(h) 
and  the  variance  is 

o'  =  1  +  phu  -  u^. 
Higher  moments  are  given  by 

n-1 


2, i/2,  ..n-l-i,h-l 


u'  =  (n-l)u'  2  +  uC  l  (1-P  )  7  '(ph) 

n  i=o 


v  i  )E(Z1 )] 


(4) 


(5) 


where  =  E( Y^| )  and  where  Ed1)  =  0  if  i  is  odd,  is  equal  to  one  if  i  is 
zero  or  two,  and  is  equal  to  (i-l)!2"^1"2^2[{(i-2)/2}!]“^  otherwise.  Ex¬ 
pressions  (3)  -  (5)  may  be  obtained  from  Johnson  and  Kotz  (1972,  pp.  112, 
114). 

Although  the  distribution  of  differs  from  the  normal,  unless  jp| 

is  near  one,  the  distribution  of  Y^  can  be  adequately  approximated  by  a 

normal  pdf.  Let  denote  a  normal  random  variable  with  mean  given  by 

(3)  and  variance  given  by  (4).  Figure  1  gives  a  comparison  of  the  pdf's 

of  Y^  and  for  h  =  0  and  p  =  .7,  9,  and  .95.  We  have  also  computed 

3/2 

the  coefficient  of  skewness  y-j  =  and  the  measure  of  kurtosis 

Y2  =  P4/U2  -  3,  where  yn  =  E[(Yh~y)n].  For  jpj  <  .5,  the  pdf's  of  Yh  and 
Z^  are  virtually  indistinguishable.  However,  as  |p|  approaches  1,  the 
pdf  of  Y^  approaches  the  pdf  of  a  normal  truncated  above  h  if  p  >  0  and 
below  -h  if  p  <0. 

[Figure  1  here] 

The  similarity  of  the  distribution  of  Y^  to  a  normal  distribution 
has  implications  for  screening  problems.  For  example,  consider  the  situa¬ 
tion  where  an  aircraft  part  is  inspected  for  fractures  using  X-rays. 
Suppose  lifetime  (L)  and  maximum  crack  length  (MCL)  are  modeled  using  a 
bivariate  normal  distribution  with  p  =  -.7.  The  average  and  standard 
deviation  (SD)  for  L  are  taken  to  be  3000  hours  and  1000  hours,  respec¬ 
tively.  Suppose  that  we  reject  any  part  with  a  crack  exceeding  1  SD  above 
average  for  MCL  (i.e.,  we  reject  about  16%  of  the  parts).  Then  using  (3) 
and  (4)  (with  h  =  1,  p  =  -.7),  y  =  .201  and  a  =  .905  and  the  distribution 
of  lifetime  for  accepted  parts  has  a  mean  of  8000  +  1000y  =  8201  hours 


and  an  SD  of  1000c  =  905  hours.  The  fifth  percentile  for  L  before  selection 
was  6355  hours.  Using  the  normal  to  approximate  the  distribution  of  L  for 
the  accepted  parts,  we  find  that,  after  screening,  only  2%  of  the  parts  should 
fail  before  6355  hrs.  [Expressions  (3)  -  (5)  are  for  computing  the  moments 
of  the  standardized  variate  Y,  conditional  on  X  <  h.  The  calculations  above 
illustrate  how  to  convert  from  y  and  a  to  the  mean  and  SD  in  original  urits.] 

APPROXIMATING  F(h,  k;  p) 

Visual  comparison  of  the  distributions  of  Yh  and  Zh  suggests  that, 
whenever  jpj  is  small,  we  may  accurately  approximate  F(h,  k;  p)  by  the 
product 

B(h,  k;  p)  =  G(h)  •  G[(k-y)/a]  (6) 

[where  y  and  a  are  as  defined  in  (3)  and  (4)],  since  Pr[Y  <  k | X  <  h]  = 

Pr[Yh  <  k]  may  be  approximated  by  Pr[Zh  <  k]  =  G[(k-y)/a]. 

As  we  investigated  the  accuracy  of  B(h,  k;  p)  for  approximating 
F(h,  k;  p)  for  various  h,  k  and  p,  we  found  that  the  error  was  less  when 
h  was  negative  and  when  jk|  <  |h|.  Therefore,  to  approximate  the  bi¬ 
variate  normal  probability  F(c,  d;  r) ,  the  error  of  approximation  is 
minimized  by  the  following  scheme: 

1.  Choose  c  and  d  so  that  jc|  >  |d|.  This  can  always  be  done,  since 
F(c,  d;  r)  =  F(d,  C;  r) . 

2a.  If  c  <  0,  set  h  =  c,  k  =  d  and  p  =  r  and  approximate 
F(c,  d;  r)  by  B(h,  k;  p) ;  or, 

2b.  If  c  >  0,  set  h  =  -c,  k  =  d  and  p  =  -r  and  approximate 
F(c,  d;  r)  =  G(d)  -  F(-c,  d;  -r)  by  G(k)  -  B(h,  k;  p). 

For  example,  to  approximate  the  probability  F(.2,  1;  .7),  set  c  =  1, 
and  following  2b,  h  =  -1,  k  =  .2  and  p  =  -.7.  Then  y  =  1.0676,  a  =  .77946, 


B(-1,  2;  -.7)  =  G(-l)  •  G[(.2-p)/a]  =  .02108,  and  the  approximation  is  G(.2)  - 
.02108  =  . 558 i 8.  The  exact  value  is  .55842,  so  our  approximation  error  is 
only  .00024. 

Using  this  scheme,  the  maximum  error  incurred  in  approximating  F(c,  d;  r) 
for  any  c  and  d 

Maximum 


r  1 

| error 

.1 

.00001 

.2 

.00006 

.3 

.00019 

.4 

.00043 

.5 

.00079 

.6 

.00121 

.7 

.00164 

.8 

.00282 

.9 

.00765 

EXAMPLE 

In  closing,  we  illustrate  the  usefulness  of  this  approximation  with 
the  example  of  Lipow  and  Eidemiller  (1964).  Rocket  motor  cases  are  known 
to  have  a  strength  S  that  is  normally  distributed  with  mean  700  pounds/ 
inch  (psi)  and  an  SD  of  100  psi.  It  is  also  known  that  peak  rocket 
operating  pressure  P  is  normally  distributed  with  mean  500  psi  and  SD 
100  psi.  (P  and  S  are  independently  distributed.)  When  a  rocket  is  fired, 
if  the  peak  pressure  exceeds  case  strength,  i.e.,  if  a  =  S-P  is  negative, 
the  case  will  rupture  and  the  rocket  will  fail.  Thus,  the  proportion  of 
cases  that  will  rupture  when  fired  is  Pr[A  <  0]  =  G(-  /2)  =  .0786  (since 
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A  is  normally  distributed  with  mean  200  psi  and  SD  141.42  psi).  In  order 
to  decrease  the  frequency  of  failure  due  to  rupture,  the  motor  cases  are 
pre-tested  to  a  pressure  of  600  psi,  and  discarded  if  they  rupture.  Hence, 
the  probability  of  rupture  for  a  case  which  passes  the  test  is  Pr[A  <  0, 

S  >  600]/Pr[S  >  600].  Now  Pr[A  <  0,  S  >  600]  =  Pr[Y  <  -/2,  X  >  -1]  with 
p  =  l/y’?,  which  also  equals  Fi-v^,  1;  -p).  To  approximate  F(-/2~,  1;  -p), 
we  compute 

B(-/2,  1;  -  .7071)  =  6(-^)6[{l  -  1.319)/. 761]  =  .0265. 

From  this  bivariate  probability,  we  can  easily  construct  the  table  of  prob¬ 
abilities: 


S  <  600 

S  >  600 

Totals 

A  <  0 

.0521 

.0265 

.0786 

A  >  0 

.1066 

.8148 

.9214 

Totals 

.1587 

.8413 

1.0000 

Hence,  the  proportion  of  failures  for  the  cases  that  passed  inspection  is 
.0265/. 8413  =  .0315.  By  screening,  we  have  reduced  the  failure  rate  from 
7.86%  to  3.15%.  However,  10.66%  of  the  motor  cases  would  have  fired 
successfully,  but  were  discarded  after  pre-testing. 

The  true  value  for  F(-v^,  1;  -.7071)  is  .0267,  so  the  error  of  the 
above  approximation  was  .0002.  This  example  illustrates  both  the  simplicity 
and  accuracy  of  the  approximation  for  bivariate  normal  probabilities  pre¬ 
sented  in  this  paper. 
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